def vv_int(a = [], dt = 1.0, vstart = 0, xstart = 0):
	"""Takes a discrete sampled acceleration vector and finds the velocity and
	displacement using Velocity Verlet method. The velocity and displacement
	vectors are assumed to have a zero initial condition by default, passing
	another value for vstart and xstart will alter the inital condition for
	velocity and displacement respectively.
	
	a: one-dimensional acceleration vector of the form a = []
	dt: time step, default = 1.0
	vstart: velocity initial condition, default value = 0
	xstart: displacement initial condition, default value = 0
	
	More info http://www.fisica.uniud.it/~ercolessi/md/md/node21.html
	
	Note: Seems to fail for -cos(x) and -sin(x), other failing cases unknown
	
	Docstring Test:
	>>> dt = 1.0
	>>> a = [9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007, 9.8000000000000007]
	>>> v, x = vv_int(a, dt, 0, 0)
	>>> print v
	[0, 9.8000000000000007, 19.600000000000001, 29.400000000000002, 39.200000000000003, 49.0, 58.799999999999997, 68.599999999999994, 78.399999999999991, 88.199999999999989, 97.999999999999986, 107.79999999999998, 117.59999999999998, 127.39999999999998, 137.19999999999999, 147.0, 156.80000000000001, 166.60000000000002, 176.40000000000003, 186.20000000000005, 196.00000000000006, 205.80000000000007, 215.60000000000008, 225.40000000000009, 235.2000000000001, 245.00000000000011, 254.80000000000013, 264.60000000000014, 274.40000000000015, 284.20000000000016, 294.00000000000017, 303.80000000000018, 313.60000000000019, 323.4000000000002, 333.20000000000022, 343.00000000000023, 352.80000000000024, 362.60000000000025, 372.40000000000026, 382.20000000000027, 392.00000000000028, 401.8000000000003, 411.60000000000031, 421.40000000000032, 431.20000000000033, 441.00000000000034, 450.80000000000035, 460.60000000000036, 470.40000000000038, 480.20000000000039, 490.0000000000004, 499.80000000000041, 509.60000000000042, 519.40000000000043, 529.20000000000039, 539.00000000000034, 548.8000000000003, 558.60000000000025, 568.4000000000002, 578.20000000000016, 588.00000000000011]
	>>> print x
	[0, 4.9000000000000004, 19.600000000000001, 44.100000000000001, 78.400000000000006, 122.5, 176.40000000000001, 240.09999999999999, 313.60000000000002, 396.90000000000003, 490.0, 592.89999999999998, 705.59999999999991, 828.09999999999991, 960.39999999999986, 1102.4999999999998, 1254.3999999999999, 1416.0999999999999, 1587.5999999999999, 1768.8999999999999, 1960.0, 2160.9000000000001, 2371.6000000000004, 2592.1000000000004, 2822.4000000000005, 3062.5000000000005, 3312.4000000000005, 3572.1000000000008, 3841.6000000000008, 4120.9000000000005, 4410.0000000000009, 4708.9000000000015, 5017.6000000000013, 5336.1000000000013, 5664.4000000000015, 6002.5000000000018, 6350.4000000000024, 6708.1000000000022, 7075.6000000000022, 7452.9000000000024, 7840.0000000000027, 8236.9000000000033, 8643.600000000004, 9060.100000000004, 9486.4000000000051, 9922.5000000000055, 10368.400000000005, 10824.100000000006, 11289.600000000006, 11764.900000000007, 12250.000000000007, 12744.900000000007, 13249.600000000008, 13764.100000000008, 14288.400000000009, 14822.500000000009, 15366.400000000009, 15920.100000000009, 16483.600000000009, 17056.900000000009, 17640.000000000007]
	"""
	
	v = [vstart]  # Velocity Initial condition with value vstart.  Default initial condition is zero (0); velocity vector
	x = [xstart]  # Displacement Initial condition with value xstart.  Default initial condition is zero (0); displacement vector
	
	for index in range(1, len(a)): # Start with one -- this allows the initial condition for v, x with "a" starting it off since it is known
		x.append(x[index-1] + (v[index-1] + .5 * a[index-1] * dt) * dt)
		v.append(v[index-1] + .5 * (a[index] + a[index-1]) * dt)
	return v, x

if __name__ == "__main__":
	from pylab import *
	import doctest
	doctest.testmod()
	
	# Use this to create an arbitrary acceleration vector, and also have an arbitrarily scaled x-axis on plots (t = [])
	dt = 1.0  # time step, made float to avoid any current python integer math problems in division
	a = []
	t = []
		
	for i in range(61):	#the range used is x - 1   #for print range(61) should show {0, 1, 2, ..., 59, 60}
		a.append(9.8)	# Function definition here, in terms of the for loop variable "i"
		t.append(i*dt)	# *dt to scale in terms of dt 

	v, x = vv_int(a, dt, vstart = 0, xstart = 0)

	# Use this to print all calculated values
	#for index in range(len(v)):
		#print "%0.2f\t%0.3f\t%0.3f\t\t%0.3f" % (index * dt, a[index], v[index], x[index])
	
	# Use this to print only the last value
	#print "%0.2f\t%0.3f\t%0.3f\t%0.3f" % ((len(a)-1) * dt, a[-1], v[-1], x[-1])
	
	# Plot the results
	subplot(221)
	title('Acceleration')
	xlabel('Time (s)')
	ylabel('Acceleration (m/s^2)')
	plot(t, a, 'b-')
	
	subplot(222)
	xlabel('Time (s)')
	ylabel('Velocity (m/s)')
	title('Velocity')
	plot(t, v, 'g-')
	
	subplot(223)
	xlabel('Time (s)')
	ylabel('Displacement (m)')
	title('Displacement')
	plot(t, x, 'r-')
	
	subplot(224)
	xlabel('Time (s)')
	title('All Plots')
	plot(t, a, 'b-', t, v, 'g-', t, x, 'r-')
	
	grid()
	subplots_adjust(left=0.11, bottom=0.07, right=.98, top=.94, wspace=.26, hspace=.28)
	savefig('verlet_out.png')
	show()